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We develop a simple, mean-field-like theory for the normal phase of a unitary Fermi gas by deriving 
a self-consistent equation for its self-energy via a momentum-dependent coupling constant for both 
attractive and repulsive universal fermions. For attractive universal fermions in the lower branch of 
a Feshbach resonance, we use zero-temperature Monte Carlo results as a starting point for one-step 
iteration in order to derive an analytical expression for the momentum-dependent self-energy. For 
repulsive universal fermions in the upper branch of a Feshbach resonance, we iteratively calculate the 
momentum-dependent self-energy via our self-consistent equation. Lastly, for the case of population 
imbalance, we propose an ansatz for higher order virial expansion coefficents. Overall, we find that 
our theory is in good agreement with currently available, high temperature experimental data. 



I. INTRODUCTION 

Fermion interactions play a central role in a wide 
range of physical systems. When these interactions are 
sufficiently strong, the physical properties of systems 
such as high-T c superconductors, neutron stars and cold 
atomic gases begin to exhibit universal behavior. Thanks 
to recent experimental progress, atoms and molecules 
can be trapped and cooled within optical or magnetic 
traps, which provide clean and highly-controllable envi- 
ronments where cold atomic systems can be studied. One 
such system, the two-component Fermi gas, is an excel- 
lent candidate for the observation of strongly interacting 
phenomena. We note that for cold atomic gases, the in- 
teraction between atoms is varied by manipulating their 
many-body bound states via the technique known as Fes- 
hbach resonance (FBR). 

At characteristic densities and ultracold temperatures, 
only isotropic and short-range s-wave scattering between 
particles can take place. This scattering can be char- 
acterized by a single parameter, the s-wave scattering 
length a. Experimentally, the s-wave scattering length 
can be tuned by using FBR [1-8]. In the unitarity limit, 
where a is tuned to ±00, the system is strongly inter- 
acting and its physical properties are independent of the 
shape of the inter-particle potential. As such, the system 
is expected to exhibit universal properties [TU] since 
its corresponding equilibrium properties depend only on 
the scaled temperature T/Tp (which is set by the energy 
scale E F and length scale I), where T F , E Fl and I are 
the Fermi temperature, the Fermi energy, and the inter- 
particle distance, respectively. Thus, by studying unitary 
Fermi gases, one learns about the equation of state for 
strongly interacting systems in general. 

Thanks to the elegant, but simple, mean-field-like the- 
ories proposed by Eagle and Legget and Nozieres and 
Schmitt-Rink, zero-temperature superfluid properties are 
well understood (TTT - FT5] . In the attractive interaction 
regime (where a < 0), atoms form Bardeen, Cooper, and 
Schreifer (BCS) pairs such that the ground state is a BCS 
superfluid. In the repulsive interaction regime (where 



a > 0), the atomic potential supports a two-body molec- 
ular bound state in vacuum such that the ground state 
is a Bose-Einstein Condensate (BEC) of these molecules. 
In between these two ground states, there is a smooth 
crossover where a changes sign as it passes through 00. 
Since the thermodynamics also evolve smoothly from the 
BCS limit to the BEC limit, a good understanding of 
the crossover regime comes from various interpolating 
schemes between these two limits [H]. However, due to 
the lack of theoretical techniques for taking into account 
strong interaction effects, the study of finite tempera- 
ture, non-superfluid (or normal) phases are challenging. 
In particular, the quantitative theoretical understanding 
of the strongly correlated problem in general is limited 
by the absence of any small physical parameters within 
the unitarity limit. While several heavily numerical ap- 
proaches [TSHin] have been used to resolve this issue, the 
simplicity of mean-field theory has been lost in the pro- 
cess. 

In this paper, we develop a simple, mean-field-like the- 
ory for strongly interacting, normal-phased Fermi gases 
at the unitarity limit. To do so, we begin in Section II 
by constructing a self-consistent theory to determine the 
self-energy of a spin-balanced, unitary Fermi gas, which 
can be accomplished by calculating the total energy, en- 
tropy and variational occupation numbers of the system. 
Continuing in Section II, we explore the thermodynam- 
ics of both the "upper branch" and the "lower branch" of 
a FBR by calculating the finite temperature equation of 
state, such that the pressure and the entropy of the sys- 
tem can be extracted and compared with experimental 
data. Finally, in Section III, we generalize an accurate 
theory concerning population balanced fermions to the 
case of population imbalanced fermions at the unitar- 
ity limit. This is done by writing the virial coefficients 
for population imbalanced fermions in terms of the virial 
coefficients for population balanced fermions, thereby al- 
lowing us to determine the grand thermodynamic po- 
tential for a system of population imbalanced, unitary 
fermions. 
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II. A SELF-CONSISTENT THEORY FOR 
UNIVERSAL THERMODYNAMICS 

In general, interaction effects in Fermi gases can be 
well described by including the appropriate self-energy 
term [H] in an expression for the total energy of the 
particles. Despite the presence of strong interaction, it 
has been argued that such a self-energy exists for uni- 
tary Fermi gases [22], however, it is impossible to derive 
an expression for it from the first principles of a Fermi 
gas alone. The main difficulty with constructing a self- 
energy for a degenerate, unitary Fermi gas is related to 
the strong interaction associated with an infinite s-wave 
scattering length. To circumvent this problem, we start 
by writing the self-energy S(fc) in terms of a distribu- 
tion function and a momentum-dependent coupling 
constant gkk' as: 



(1) 



where V is the system volume. Here, the distribution 
function is taken as a variational function. In this 
Hartree-Fock-like self-energy, gkk' is given by: 
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where the phase shift 5 and the scattering length a are re- 
lated to one another by: 6 = — arctan(|& — k'\a/2). Here, 
m is the mass of a particle and \ k— k'\/2 = q is the relative 
momentum between a particle and the scattering parti- 
cle. Note that in the limit where (k — k')\a\/2 <C 1, the 
momentum-dependent coupling constant transforms into 
the ordinary momentum-independent, mean- field cou- 
pling contasnt g — ATrh 2 a/m, while at the unitarity limit 
(where a — > ±oo), <5 — ^ <So = Tt/2, a constant value. 
Now, by calculating the total energy: 
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the entropy: 



S = -^[n fc lnn fc + (1 - n fe )ln(l - rife)] 



(4) 



and the number of particles: N — 2^2 k nk, the grand 
thermodynamic potential Q = E — TS — fiN is derived. 
Then, by using the relation: N = —dil/d/i, the varia- 
tional occupation number nk is obtained as: 



(<*-/*) -f 1 
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where efe = h 2 k/2m + T,(k) is the total single particle en- 
ergy of a particle with momentum k and fi is the chemical 
potential of the system. Lastly, by inserting this occu- 
pation number and the unitary value of the momentum- 
dependent coupling constant into Eq. ([I]), we complete 
the derivation of a self-consistent equation for the self- 
energy £(fc) at unitarity: 
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However, since real-space densities for ultracold Fermi 
gases are so dilute, they also have large momentum-space 
densities associated with them. As a result, we can con- 
vert the sum in Eq. ^ into an equivalent integral over all 
of momentum-space. Then, by defining the dimcnsionless 
variables: = -q, /3£(fc) = V{^y) and 7 = Hk\J /3/(2m) 
(where j3 = l/fc^T), we have: 



V( 7 ) = 
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By expanding — VI m terms of spherical harmonics 
and using their orthonormal properties, we can solve the 
angular part of Eq. ^ for both the 7 < 7' case and the 
7 > 7' case. Hence, we obtain: 
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for 7 < 7' 
for 7 > 7' 
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Finally, we distribute Eq. ^ into two, separate integrals 
(corresponding to the two cases) and simplify the results 
to arrive at the concluding form of our self-consistent 
equation for the self-energy V(pf) at unitarity: 
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Upper Branch Thermodynamics of a Feshbach 
Resonance 

The BEC side of a Feshbach resonance is characterized 
by a positive scattering length (a > 0) and by a ground 
state which involves bound pairs of fermionic molecules in 
condensate. This is known as the "lower branch" of a Fes- 
hbach resonance. By contrast, in the "upper branch" of 
a Feshbach resonance, the system's wavefunction consists 
of scattering states such that we can neglect these bound 
pairs and their corresponding binding energies |23j . Mo- 
tivated by a series of recent experiments conducted in this 
(metastable) upper branch state [2~il - f2l)] . static and dy- 
namic properties of Fermi gases have been theoretically 



FIG. 1: The self-energy of a unitary Fermi gas in the 
upper branch of a Feshbach resonance as related to its 
momentum. Values of 7 = 0.1, 1 and 2 are shown in 
black, dark gray and light gray lines, respectively. 

studied [27l431| . To further explore the thermodynamics 
of the upper branch, we set 5$ = —n/2 (i.e. its limiting 
value as a — > +00) in order to find the numerical solu- 
tion of Eq. ([9]) iteratively. We note that while Eq. ^ 
converges very rapidly at first, it requires more iterations 
at relatively small values of 7. As such, the calculated 
self-energy V (7) for the upper branch is shown in FIG. [I] 
for three different values of 7, while the calculated oc- 
cupation numbers 71(7) for the upper branch are simi- 
larly shown in FIG. [2j Additionally, we show the calcu- 
lated pressure P* = P(£)A 3 /fcsT for the upper branch in 
FIG. [3] where £ = exp(— /i/feT). It is important to note 
that the tail of the momentum distribution should decay 
with 7~ 4 behavior as related to the Tan Contact [3"2"H3"1] . 
However, the Contact density as a function of temper- 
ature was recently determined by Boettcher et al. for 
the self-energy of ultracold fermions in the BEC-BCS 
crossover by using methods in non-perturbative quantum 
field theory [35] . Their findings indicate that the Contact 
density is not a monotonic function of temperature, and 
that its maximum occurs at approximately 1.25T C (where 
T c is the critical temperature for phase transition). Ad- 
ditionally, their results show that the Contact density 
becomes very small for temperatures greater than 2T C , 
which is reinforced by the findings of Enss and Hauss- 
mann, who have determined that the Contact density for 
a unitary Fermi gas is C = 0.086A;|. at T = Q.5T F [37)] . 
Since our theory is valid for the normal phase of a Fermi 
gas at unitarity, we believe that the Contact plays a very 
small role in the actual decay of our momentum distri- 
bution function, especially at higher temperatures. 



BCS Side Thermodynamics of a Feshbach Resonance 

The BCS side of a Feshbach resonance is characterized 
by a negative scattering length (a < 0) and by a ground 



FIG. 2: The occupation numbers of a unitary Fermi gas 
in the upper branch of a Feshbach resonance as related 
to momentum. Values of 7 = 0.1, 1 and 2 are shown in 

black, dark gray and light gray lines, respectively. 
(Dashed lines of the same hue are the non-interacting 
occupation numbers corresponding to these 7 values.) 
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FIG. 3: The pressure of a unitary Fermi g as as a 
function of £ in the upper branch of a Feshbach 
resonance. Here, P* = P(£)A 3 /ksT, where 
£ = exp(-^/k B T), 

state which involves Cooper pairs of fcrmionic molecules 
in condensate. To further explore the thermodyanmics of 
this lower branch state, we set S = +ir/2 (i.e. its limiting 
value as a — > —00) in order to find the numerical solution 
of Eq. ([9| iteratively. However, upon doing so, we realize 
that Eq. ([9]) does not converge for relatively small values 
of 7 as it did before (for small values of £). In principle, 
we could introduce a lower cutoff value for the momentum 
to get around this problem, but our resulting numerical 
iteration may not be reliable. As such, we choose to solve 
Eq. ^ via one-step iteration instead. The motivation 
for doing so is two-fold: not only do we want to get away 
with the non-convergency at low values of momenta, but 
we also want to find an approximate analytical expres- 
sion for the self-energy. However, obtaining an accurate 
expresssion for the self-energy requires us to first choose 
an accurate starting point for our one-step iteration (i.e. 
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the zeroth order self-energy). 

In general, the unitarity Fermi gas at zero-temperature 
has been studied with both the heavily numerical Monte 
Carlo method |18j and with renormalization group the- 
ory [37]. Monte Carlo calculations at zero-temperature 
have shown that the self-energy is given by: fiE(k) = Afi 
(where A ~ —0.4045 16, 38] El]), while renormalization 
group theory has shown that the self-energy has a weak 
momentum and frequency dependence at the unitarity 
limit [37]. Based on these zero-temperature theoretical 
results, we make the ansatz ?j£(fc) = Afi as the zeroth 
order self energy. Thus, for one-step iteration, we define: 
h = rf — V = i] — Arj. Doing so (while expanding the 
denominator) allows us to write the first integral on the 
right hand side of Eq. ([£]) as: 



n=0 
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Completing these two integrations yeilds the expression: 
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where g n {x) 
and f n (x) = 

The two functions: Erf{x) and Erfi(x), are the usual 
error and imaginary error functions, respectively. Note 



that for h < 0, only the second term in Eq. (Ill con- 



tributes, whereas for h > 7, only the first term con- 
tributes. Thus, for h < 0, we must set f n {h —¥ 0) =0 in 
the second term. The second integral on the right hand 
side of Eq. Q can be evaluated explicitly to obtain: 



h = (ln[e' 1 + e? ] - 7 2 )/2 



(12) 



Finally, within one-step iteration, we find an analyti- 
cal expression for the momentum-dependent, finite tem- 
perature self-energy as: V = — 2(^/7 + I 2 ). Note that 
while we have completed the momentum integration, the 
self-energy now has the form of an infinite converging 
series. Now, in order to test the accuracy of our itera- 
tive analytical theory, we calculate the finite temperature 
equation of state and then compare our results with ex- 
perimental data. In the unitary limit, the pressure can be 
written in a universal form as: P(p,T) = Pi(/i, T)/i(£). 
Here, Pi(/i, T) is the pressure of a single-component, non- 
interacting Fermi gas: 



k B T 2 
A 3 V^Jo 



Viln[l + z^e'^dt (13) 




FIG. 4: The pressure of a unitary Fermi g as as a 
function of £ in the lower branch of a Feshbach 

resonance. Here, h(£) = P{pL,T)/P^(jjb,T), where 
£ = exp(— h/UbT). Our theory (black line) vs. 
experimental data |3jJ| (gray points). 



where A = y/2nh 2 /mksT is the thermal wavelength and 
£ = cxp(~n/kBT). Upon completion of the self-energy 
calculation, we compute the system pressure via the rela- 
tion: P(n, T) = —dfl/dV, and then extract the universal 
function /i(£). This function is plotted in FIG. [4] together 
with experimental data from Nascimbene et al. [40] . As 
can be seen, our theory is in reasonable agreement with 
experimental values at higher temperatures, however, we 
begin to see a noticeable deviation at very low tempera- 
tures. This is due to the fact that our theory is only valid 
for the normal phase, wheras the experimental system is 
in the superfluid phase. 

We note that while the jump in Sq by tt at unitarity 
will correspond to a jump in our self-energy, the thermo- 
dynamic quantities of the system will remain continuous 
throughout the BEC-BCS evolution. This is due to our 
treatment of the upper branch on the BEC side, where we 
have chosen to neglect the binding energy of ground-state 
pairs. As such, the jump in our self-energy should vanish 
if we were to take these binding energies into account. 
However, since this was not the case, FIG. [3] and FIG. [4] 
are qualitatively the same, but quantitatively different. 
Note that this behavior has already been experimentally 
verified via upper branch energy measurements |24) , thus 
we believe that these two figures would also be quanti- 
tatively the same if the jump in our self-energy was not 
present. Regardless, we have chosen to validate our the- 
ory by using data [40j from experiments performed on the 
BCS side of a Feshbach resonance, where its numerical 
predictions are most suited for comparison with thermo- 
dynamic measurements. 

Similarly, the temperature dependence of the entropy 
and energy for harmonically trapped, unitary fermions 
is measured by Thomas's group at Duke University [41] . 
To investigate this, we used the local density approxima- 
tion (LDA) to evaluate the entropy and energy of these 
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trapped fermions. In LDA, the local chemical potential 
(i is written in terms of the central chemical potential 
Ho and the trapping potential V{r) = m[u)\{x 2 + y 2 ) + 
to 2 z 2 ]/2 as: \x = /i - muj 2 r 2 /2, where ui — (a^u^) 1 / 3 
is the average trapping frequency. In doing so, the total 
number of particles: 

N = [ d 3 rn(r) = — ^ f P(r)dr (14) 
J muj 2 J 

and the total energy: 

E = 12ir J r 2 P(r)dr (15) 

can easily be converted into integrals over the chemical 
potential as: 




1 2 3 4 5 6 

S/Nk e 

FIG. 5: The entropy as a function of energy for 
harmonically trapped, unitary fermions. Our theory 
(black line) vs. experimental data [H] (gray points). 



N 



E = 
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(16) 
(17) 



where r\ = fiji and 770 = /3/Uq. The entropy is then given 
by: S = 4E(3/3 - r) N. Now, by defining the Fermi 
energy and Fermi temperature of an ideal Ferm i ga s as: 
Ep = ksTp = (SNy^Huj, we can combine Eq. (16) and 



(17) to yeild the expression: 
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Hence, for given values of 770 and T, we can solve the 
above set of equations for the entropy and the energy. 
As such, the calculated entropy as a function of energy 
is shown in FIG. [5] along with experimental data taken 
from Luo and Thomas 1411. 



virial coefficients are temperature independent. In this 
section, we generalize the BvDM approach to the case of 
population imbalanced, unitary fermions. 

First, we summarize the original BvDM approach |46j 
by noting that the grand thermodynamic potential of a 
population balanced Fermi system can be written as: 



00 



-k B TZi(f3) 2_^(Abi)z 
1=0 



(19) 



where f^ ) is the grand thermodynamic potential of an 
ideal Fermi gas, Z\(J$) is the one-particle partition func- 
tion and Abi — bi — b[ ^ is the /,-particle cluster integral 
relative to an ideal Fermi gas. As FBR is related to the 
forming and dissolving of two-body pairs, BvDM has pro- 
posed that higher order cluster integrals are expressible 
in terms of the two-body cluster A&2- Assuming the l- 
body cluster is one particle interacting with I — 1 paired 
particles, the Z-particle cluster integral is given by: 



III. POPULATION IMBALANCED FERMIONS 
AT UNITARITY 

Recent experiments concerning population imbalanced 
fermions [32J 03] has triggered a new direction in theo- 
retical research devoted to the study of unitary fermions 
in the presence of population imbalance [32 HH |4"S] . 
For population balanced, two-component fermions at 
unitarity, R. K. Bhaduri, W. van Dijk and M. V. N. 
Murthy (BvDM) have proposed a parameter-free, high- 
temperature equation of state based on a virial cluster 
expansion |46) that shows excellent agreement with ex- 
perimental results over a wide range of fugacity. Their 
basic assumption is that higher order cluster integrals can 
be written in terms of two-particle clusters. This is jus- 
tified because only two-body scattering effects are domi- 
nant for dilute atomic gases, even at unitarity where the 



A6 i = (- 1 ) , 25r ( 2 °) 

for I > 2, where ai = (I — — 2)/2. As was mentioned 
earlier, the BvDM ansatz for a population balanced sys- 
tem of fermions shows excellent agreement with experi- 
mental results over a wide range of fugacity. Hence, for 
the remainder of this paper, we generalize this ansatz to 
the case of population imbalanced, unitary fermions. 

The grand thermodynamic potential of a population 
imbalanced Fermi system can be written as: 

00 00 

n = -k B Tz 1 (p)J2zZ b ^ z r kz i ( 21 ) 

n=l k=0 

where b n< k is the n-th virial coefficient for a configura- 
tion of n — k spin- up fermions and k spin-down fermions. 



6 



We also note that b ni k has the properties: b n , n -k = b n ,k 
and Y^k=o bn,k = b n . Now, by defining the virial coef- 
ficient difference relative to non-interacting fermions as: 



A6„, fc = b n „ 



,(o) 
Vfc> 



we have: 



n - r» = -k B TZ x {p) E( A6 ».*K~ 



kk 



(22) 



n=l fc=0 



where J7o — ^o<t« Here, we note that Oo CT is the grand 
thermodynamic potential for a non-interacting fermions, 
and is given by the expression: 




n 
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k B T 2 



]dt. 



(23) 



Since the interaction occurs only between fermions with 
opposing spins, the virial coefficients for population im- 
balanced fermions can be written in terms of virial coef- 
ficients for population balanced fermions. Furthermore, 
A^n.o = 0, and we find that: A„ fc = A n /(n — 1) for 
n > 2. Finally, by putting everything together, the grand 
thermodynamic potential for a population imbalanced, 
unitary Fermi system can be written as: 



n-n = -k B Tzm £ U n _u E *r k 4 ( 24 ) 

n=2 y I k=l 

Hence, by using the relation: P = —dtl/dV, we see that 
the pressure for a population imbalanced, unitary Fermi 
system can be similarly written as: 



k B T ~ (-1)»A&3 

n=2 



■(n-1) 
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Lastly, we extract the universal function £) from Eq. 



( 25 1 to obtain: 



h( v ,0 = i 
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where r/ = fii/fi^ is the ratio between the two chemical 
potentials. We note that when rj = 1.0, the chemical po- 
tentials are equal and Eq. ( 26 ) reduces to its population 
balanced form as depicted in FIG. [4] By contrast, h(r],^) 
is plotted in FIG.[6]as a function of for three different 
values of rj, along with experimental data 40J pertaining 
to the population balanced case for verification purposes. 



IV. CONCLUSIONS AND REMARKS 

Two important points of interest to consider are the re- 
cent experimental findings by Zwierleins's group |47j and 



FIG. 6: The finite temperature equation of state for a 
population imbalanced, unitary Fermi gas. Values of 
r) = 1.0, 0.5 and 0.25 are shown in black, dark gray and 
light gray lines, respectively. The gray points represent 
experimental data [40] for the population balanced case, 
which is recovered from Eq. (261 when 77 = 1.0. 



the recent theoretical findings by Van Houcke et al |48j . 
Zwierlein's group observed the superfluid phase transi- 
tion in a strongly-interacting Fermi gas by using high- 
precision measurements of the local compressibility, den- 
sity and pressure. Their data was able to completely de- 
scribe the universal thermodyanmics of such Fermi gases 
without the use of any fit or external thermometer. Sim- 
ilarly, Van Houcke et al. computed and measured the 
equation of state for a normal, unitary Fermi gas. Their 
data showed excellent agreement with their theory that 
a series of Feynman diagrams can be controllably re- 
summed in a non-perturbative regime using a Bold Dia- 
grammic Monte Carlo approach. Wc note that while the 
newer data from Zwierlein's group is highly accurate [4 7) . 
its normal phase measurements are nearly identical to the 
measurements made by Nascimbene et al. [4"0] . 

In conclusion, we have presented a self-consistent the- 
ory to determine the self-energy of a strongly interact- 
ing, normal-phased Fermi gas at unitarity. We have 
also shown that this self-energy can be used to calcu- 
late universal thermodynamic properties of a Fermi gas 
for both the upper branch of a Feshbach resonance and 
the lower branch of a Feshbach resonance. In addition, 
we have demonstrated that higher order virial expan- 
sion coefficients for population imbalanced fermions can 
be written in terms of virial expansion coefficients for 
population balanced fermions, which makes calculating 
the grand thermodynamic potential for population im- 
balanced Fermi systems at unitarity a much less cumber- 
some task. Overall, we find that our theory is in good 
agreement with currently available experimental data, 
indicating that there may be promising advancements 
ahead for further theoretical research regarding fermion 
interactions at the unitarity limit. 
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